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OF  THE  EULER  PROBLEM  FOR  SWEPT  WINGS* 

by 
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/  SUMMARY 

h  - 

A  multigrid  cell-vertex  finite  volume  Euler  method  has  been 
used  to  calculate  steady  inviscid  transonic  flow  past  the  ONERA  M6 
wing.  The  treatment  of  the  far-field  boundary  conditions  includes 
the  effect  of  velocity  perturbations  generated  from  Klunker's 
analytic  asymptotic  solution  to  the  transonic  small-disturbance 
equation.  The  geometry  of  the  wing  tip  is  modelled  with  three  suc¬ 
cessively  finer  C-0  grids.  The  results  obtained  show  marked  differ¬ 
ences  in  comparison  with  those  obtained  on  C-H  grids.  They  indicate 
that  a  shock-wave  mechanism  can  contribute  to  the  generation  of  tip 
vortices  in  compressible  inviscid  flow,  and  that  for  reliable  pre¬ 
diction  of  wing  performance  it  is  important  to  model  well  both  the 
wing-tip  geometry  and  the  flow  around  the  tip.  K,  ./<■  ,*  /  ■ 
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1  INTRODUCTION 

1  1 

The  work  of  Ni  has  been  extended  in  two  and  three  dimensions  by  Hall  and 

3 

Salmond  .  Numerical  solutions  to  the  Euler  equations  for  transonic  lifting  flows 
past  aerofoils  and  wings  are  obtained  by  approximating  the  steady  flux  balance 
for  each  computational  cell  using  values  of  the  flow  variables  stored  at  the 
vertices  of  each  cell.  A  Lax-Wendroff  time-stepping  algorithm  is  used  to  update 
each  variable,  with  boundary  conditions  where  appropriate,  and  these  changes  in 
the  solution  are  then  used  on  coarser  meshes  to  accelerate  the  convergence  of  the 
solution  to  the  steady  state.  This  paper  describes  some  recent  developments, 
mainly  in  the  imposition  of  far-field  boundary  conditions  and  in  the  modelling  of 
the  flow  around  the  wing  tip. 

The  usual  practice  in  calculating  the  three-dimensional  flow  past  a  wing 
has  been  to  set  conditions  on  the  far-field  boundary  equal  to  the  free-stream 
conditions.  This  assumption  is  loosely  justified  by  the  knowledge  that,  except 
in  a  wake,  a  disturbance  would  decay  more  rapidly  with  distance  in  three  dimen¬ 
sions  than  in  two.  Here,  instead,  we  replace  the  free-stream  conditions  by  a 

three-dimensional  far-field  solution.  For  'inflow'  the  conditions  are  derived 

4 

from  the  transonic  small-disturbance  approximation  of  Klunker  ,  while  for  'out¬ 
flow'  a  perturbation  pressure  related  to  the  cross-flow  in  the  numerical  solution 
is  prescribed.  The  results  are  expected  to  be  more  accurate  than  those  given  by 
assuming  free-stream  conditions  and  to  provide  an  indication  of  the  range  of 
validity  of  the  latter  assumption.  Comparisons  are  shown  for  a  range  of  distances 
from  the  wing  to  the  far-field  boundary. 

It  is  becoming  recognised  (for  example,  see  Ref  5)  that  the  resolution  of 
flow  detail  at  the  wing  tip  requires  both  the  representation  of  the  tip  geometry 
and  Che  use  of  an  appropriate  grid  around  the  tip.  We  outline  here  the  steps 
taken  in  generating  C-0  grids  from  C-H  grids  formed  by  spanwise  stacking  of  two- 
dimensional  planar  C-grids.  Results  are  shown  for  calculations  performed  on  the 
ONERA  M6  wing  using  both  C-H  grids  and  three  successively  finer  C-0  grids.  We 
begin  however  with  an  outline  of  the  basic  numerical  method. 

2  CELL  VERTEX  FINITE  VOLUME  EULER  METHOD 

Since  only  the  steady  state  is  of  interest  the  unsteady  Euler  equations  in 
three  dimensions  are  expressed  in  conservation  form  as  the  isenthalpic  system 

u^  +  fvu)^  ♦  g(u) ^  +  h  (u)  ^  *  0  , 


(!) 
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where  the  state  vector  a  and  the  flux  vectors  f(u),  §(u),  h(u)  are  given 
as  follows: 
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Here  p  is  the  density,  u,  v,  w  are  Cartesian  components  of  velocity  and  p 
is  the  pressure  given  by  Bernoulli's  equation 

P  "  ~  (I  “  -  D(u2  +  v2  +  w2))  , 

where  y  is  the  ratio  of  specific  heats,  taken  to  be  1.4.  The  equations  have 
been  non-dimensionalised  with  respect  to  a  typical  length,  stagnation  sound 
speed  and  stagnation  pressure. 

Integrating  equation  (1)  over  a  control  element  of  volume  V  and  using 
Che  divergence  theorem  yields  for  the  steady  state  the  boundary  integral  (or 
residual) 


7  =  0  ‘ 


(2) 


Dividing  our  region  of  physical  space  surrounding  the  wing  into  hexahedral  cells 
we  seek  co  approximate  equation  (2)  for  each  cell.  The  fluxes  through  each  of 
the  faces  of  a  particular  cell  are  calculated  from  the  values  of  the  flow  vari¬ 
ables  stored  at  the  vertices  of  that  cell  (see  Fig  1)  and  are  summed  to  give  the 
net  flux  into  the  cell. 


Given  an  initial  flow  the  steady  solution  u  satisfying  equation  (2)  for 
all  cells  is  obtained  by  a  Lax-Wendroff  time-marching  algorithm.  From 
equation  (I)  the  change  in  the  solution,  5un  =  un+'  -  u11  ,  in  time-step  it  is 

given  by 
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the  two  parts  of  the  change  5u  at  the  point  P  are  calculated  as  follows  (see 
Fig  2).  The  first  term  in  expression  (3)  is  given  by  averaging  the  residuals  (2) 
for  the  eight  cells  surrounding  P  .  The  second  term  is  recast  by  integrating 
over  the  inner  hexahedron  surrounding  P  formed  by  the  centres  of  the  eight 
cells.  The  divergence  theorem  again  yields  an  integral  over  the  surface  of  this 
hexahedron,  evaluated  as  before,  the  values  of  the  flow  variables  used  being 
obtained  as  averages  of  the  values  at  the  vertices  of  the  eight  cells. 


To  prevent  odd-even  point  decoupling  and  to  capture  shock-waves  an  arti¬ 
ficial  viscosity  is  added  to  5u  .  Convergence  acceleration  is  achieved  using 
local  time-stepping  and  a  multigrid  procedure.  Changes  in  the  solution  on  the 
coarse  grids  are  interpolated  linearly  to  the  finer  grids  as  previously  but  are 
then  smoothed  by  an  averaging  operator. 

3  FAR-FIELD  BOUNDARY  CONDITIONS 


In  the  present  method  the  inflow  condition  is  that  the  two  components  of 
velocity  in  the  plane  locally  tangential  to  the  boundary  are  specified  together 
with  an  isentropy  condition  for  the  density.  To  determine  corrections  to  the 
free-stream  condition  in  the  far-field  we  firstly  consider  Klunker's  expression 
for  the  transonic  small-disturbance  potential  4  at  a  point  P(x,y,z)  .  In  the 
far-field  Klunker  shows  that  the  dominant  term  is  that  representing  the  lifting 
effects  of  the  wing,  namely  the  following  integral  over  the  wing  surface  (see 
Fig  3) 
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where  C,  n  are  local  coordinates  on  the  wing  and 
,2  2  .  .  .2  ,  .,2 


2  2  2  2 
R  =  (X  -  O  +  S  (y-  c) 


3~z  with  6”  «  1  -  M  .  Au  is  interpreted  here  as  the  difference  in  pressure 
between  the  upper  and  lower  surfaces.  Due  to  deficiencies  in  the  potential  model 
4  is  only  a  good  approximation  away  from  the  wake,  and  in  fact  is  only  used  to 
generate  perturbations  to  the  velocities  given  at  'inflow'. 


4  is  differentiated  to  give  expressions  for  the  perturbations  (Su,6v,Sw) 
P 

from  the  free-stream  velocities  (u  ,v  .w^) .  The  resulting  integrals  are  evalu¬ 
ated  numerically  for  each  point  P  on  the  'inflow'  part  of  the  far-field  bound¬ 
ary.  The  modified  velocities  (u  +  Su,  v^  +  6v,  +  <5w)  are  then  used  to  cal¬ 

culate  the  two  imposed  components  of  the  tangential  velocity.  If  the  perturba¬ 
tions  are  set  to  zero  we  regain  the  free-stream  boundary  conditions. 
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At  outflow  a  boundary  condition  on  pressure  p  is  assumed  which  includes 
the  effect  of  the  component  of  velocity  q^  normal  to  the  free-stream  direction 
It  can  be  shown  that  if  variations  in  the  function  p/oY  are  small,  and  if  at 
outflow 


2  2  2 
U  +  V  +  w 


where  is  the  free-stream  velocity,  then  the  Bernoulli  equation  can  be 

approximated  by 


p 


Y/(Y-1> 


This  is  the  assumed  condition. 

4  WING  TIP  MODELLING 

The  stages  in  producing  a  C-0  grid  from  the  C-H  grid  comprising  two- 
dimensional  planar  C-grids  stacked  spanwise  are  as  follows: 

(i)  discard  the  C-grids  lying  beyond  the  wing  tip; 

(ii)  complete  the  surface  closure  for  the  M6  wing  using  circular  arcs  as 
continuations  of  the  wing  generators.  The  (finest)  surface  grid  has 
each  arc  subdivided  into  32  equal  segments; 

(iii)  as  sketched  in  Fig  4,  rotate  the  C-grid  corresponding  to  the  tip 
station  180°  about  an  axis  running  chordwise  through  the  leading  and 
trailing  edges.  Equal  subdivision  of  the  arcs  swept  out  provide  the 
field  grid.  Lines  of  semi-polar  singularity  are  left  along  the  axis 
of  rotation  forward  of  the  leading  edge  and  aft  of  the  trailing  edge. 

After  these  steps  adjustments  were  made  to  produce  the  final  grid.  These  inclu¬ 
ded  bending  and  packing  of  the  C-grids  to  ensure  smooth  spanwise  distributions, 
the  addition  of  a  further  16  C-grids  on  the  main  wing,  and  contracting  the  grid 
around  the  trailing  edge  of  the  tip  where  the  tip  radius  is  small.  Thus  from  an 
original  C-H  grid  of  256  x  32  *  48  cells  was  produced  a  new  C-0  grid  of 
256  x  32  x  64  cells.  Details  of  the  medium  grid  (128  x  16  x  32)  on  the  wing  sur 
face  near  Che  tip  are  shown  in  Figs  5  and  6. 

5  RESULTS 

The  numerical  method  was  used  in  two  forms,  appropriate  to  the  C-H  and 
C-0  grids,  to  calculate  the  flow  past  the  ONERA  M6  wing,  with  *  0.84  and 


a  ■  3.06°.  To  evaluate  the  effect  of  the  new  boundary  conditions  the  tar-field 
boundary  was  set  at  varying  distances  from  the  wing  -  1.01,  1.69  and  3.37  spans 
(3,  5  and  10  root  chords  respectively)  -  and  the  calculations  performed  with  and 
without  the  perturbations  modifying  the  free-stream.  The  effect  on  is  shown 

in  Fig  7  for  the  calculations  on  the  medium  C-H  and  C-0  grids.  With  the  free- 
stream  conditions  the  variation  in  is  around  32  while  for  the  perturbed  con¬ 

ditions  it  is  around  0.72.  All  subsequent  calculations  were  performed  with  the 
perturbed  conditions  with  the  boundaries  set  at  1.69  spans. 

Next  we  detail  the  comparison  between  the  tip  flows  on  the  medium  C-H  and 
C-0  grids.  The  spanwise  lift  distributions  in  Fig  8  show  a  marked  spike  at  the 
tip  for  the  calculation  on  the  C-H  grid,  indicating  a  vortex.  This  is  absent  for 
the  calculation  on  the  C-0  grid  and  examination  of  the  cross-flow  shows  little 
sign  of  a  vortex,  although  the  cross-flows  are  much  larger  than  for  the  C-H  grid. 
Shown  in  Figs  9  and  10  are  the  upper  surface  Mach  number  contours  in  the  tip 
region  for  the  two  cases.  As  might  be  expected  the  flows  are  very  similar  inboard 
but  quite  different  at  the  tip.  For  the  calculation  on  the  C-0  grid  the  contours 
near  the  tip  trailing  edge  indicate  the  presence  of  high  velocities  there. 

Finally  we  compare  the  calculations  performed  on  the  successively  finer 

-2/3 

C-0  grids.  Fig  11  shows  the  value  of  C^  plotted  against  (number  of  cells)  , 
the  straight  line  indicating  second-order  accuracy.  The  spanwise  lift  distribu¬ 
tions  are  plotted  in  Fig  12.  Note  the  small  kink  at  the  tip  for  the  finest  grid; 
indeed  (Fig  13)  a  plot  of  the  cross-flow  velocities  obtained  on  this  grid,  very 
close  to  the  trailing  edge  (98.52  chord),  shows  a  region  of  recirculating  flow. 

The  production  of  the  recirculating  flow,  or  tip  vortex,  can  be  explained 
by  examining  the  cross-flow  upstream.  Having  become  sonic  at  about  892  chord  the 
cross-flow  Mach  number  reaches  3.0  at  96.52  chord.  A  shock  in  the  cross-flow 
plane  can  be  seen  in  the  velocity  vectors  and  in  the  contours  of  cross-flow  Mach 
number  in  Figs  19  and  15  obtained  on  the  finest  grid.  The  corresponding  plot  of 
total  vorticity  contours  in  Fig  16  shows  that  while  some  vorticity  is  produced  at 
the  junction  of  the  circular  arc  with  the  wing  generator,  more  is  produced  through 
the  shock-wave.  This  vorticity  is  convected  downstream  and  may  therefore  be  the 
origin  of  the  vortex  shown.  Further  work  needs  to  be  done  to  assess  whether 
separation  is  occurring  in  addition. 

The  marked  differences  between  the  results  on  the  C-0  grid  and  those  on  the 
C-H  grid  have  important  practical  implications.  They  show  that  the  flow  pattern 
around  the  wing  tip  and  the  associated  spanwise  lift  distribution  are  strongly 
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misrepresented  when  calculations  are  performed  on  the  commonly  adopted  C-H 
grid.  For  reliable  prediction  of  wing  performance  the  wing-tip  geometry  and  the 
flow  around  the  tip  must  be  well  modelled,  at  least  for  wings  of  low  aspect  ratio. 

6  CONCLUSIONS 

Improved  far-field  boundary  conditions  and  wing  tip  modelling  have  been 
proposed  and  tested.  The  boundary  conditions  allow  the  distance  to  the  far-field 
boundary  to  be  reduced  considerably  with  little  change  in  lift.  Calculations  on 
C-H  and  C-0  grids  show  very  marked  differences  in  the  flow  at  the  wing  tip,  there 
being  a  distinct  tip  vortex  on  the  C-H  grids.  Although  cross-flows  are  much 
larger  on  the  C-0  grids  there  is  little  sign  of  a  tip  vortex  except  for  the 
finest  grid.  In  this  case  a  vortex  appears  to  originate  at  least  in  part  in 
vorticity  generated  by  a  shock  in  the  cross-flow.  For  reliable  prediction  of 
wing  performance  it  is  important  to  model  well  both  the  wing-tip  geometry  and 
the  flow  around  the  tip. 

Acknowledgments 

The  authors  thank  Cray  Research  (UK)  Ltd  and  AERE  Harwell  for  providing 
assistance  and  computing  facilities  essential  for  the  completion  of  the  work 
described. 


TM  Ac  21 }9 


REFERENCES 


No,  Author  Title,  etc 

1  R.-H.  Ni  A  multiple  grid  scheme  for  solving  the  Euler  equations. 

AIAA  «7.,  20,  11,  1565-1571  (1982) 

2  M.G.  Hall  Cell-vertex  schemes  for  solution  of  the  Euler  equations. 

Proc.  Conf.  on  Num  Meth  for  FI  Dyn,  University  of 

Reading,  1985,  (Eds.  K.W.  Morton  and  M.J.  Baines)  Oxford 
University  Press,  303-345  (1986) 

3  D.J.  Salmond  A  cell-vertex  multigrid  scheme  for  solution  of  the  Euler 

equations  for  transonic  flow  past  a  wing. 

Proc.  ICNMFD  10,  Beijing,  China,  1986  (Eds.  F.G.  Zhuang 
and  Y.L.  Zhu)  Springer-Ver lag  Lecture  Notes  in  Physics 
No. 264,  549-553  (1986) 

4  E.B.  Klunker  Contribution  to  methods  for  calculating  the  flow  about 

thin  lifting  wings  at  transonic  speeds  -  analytic 
expressions  for  the  far-field. 

NASA  TN  D-6530  (1971) 

5  P.  Sacher  Numerical  solutions  for  three-dimensional  cases  -  swept 

wings. 

Chapter  7,  "Test  cases  for  inviscid  flow  field  methods". 
AGARD  Advisory  Report  No. 211  (1985) 


Figs  5-10 


Fig  5  Surface  grid  at  the  wing  tip 
(leading  edge) 


Fig  6  Surface  grid  at  the  wing  tip 
(trailing  edge) 
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Fig  11  Variation  of  Cl  with  refinement 
of  the  C-0  grid 
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Fig  13  Cross-flow  velocity  vectors  at 

the  wing  tip  very  near  the  trailing 
edge  (98.5%  chord) 
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Fig  t5  Cross-llow  Mach  contours  at 
the  wing  tip  near  the  trailing 
edge  (96.5%  chord) 
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Fig  12  Variation  of  the  spanwise  lift 
distribution  with  refinement 
of  the  C-0  grid 


Fig  14  Cross-flow  velocity  vectors  at 
the  wing  tip  near  the  trailing 
edge  (96.5%  chord) 


Fig  16  Vorticity  (|fi|)  contours  at  the 
wing  tip  near  the  trailing  edge 
(96.51S  chord) 
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A  multigrid  cell-vertex  finite  volume  Euler  method  has  been  used  to  calculate 
steady  inviscld  transonic  flow  past  the  ONERA  M6  wing.  The  treatment  of  the  far- 
fleld  boundary  conditions  includts  ths  sffect  of  velocity  perturbations  generatsrd 
fror  Klunker's  analytic  asymptotic  solution  to  the  transonic  smalWlsturbancs 
equation.  The  geometry  of  the  wing  tip  is  modelled  with  three  successively  finer 
C-0  grids.  The  results  obtained  show  marked  differences  in  comparison  with  those 
obtained  on  C-ll  grids.  They  indicate  that  a  shock-wave  mechanism  can  contribute 
to  the  generation  of  tip  vortices  in  compressible  inviscld  flow*  and  that  for 
rf liable  prediction  of  wing  performance  it  is  important  to  model  well  both  the 
wing-tip  geometry  and  the  flow  around  the  tip. 
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